function productivitydraws = Pareto_draw(beta,theta,qual_vec,draws)
%Make a function for drawing from a Pareto distribution

CDFvals = rand(draws,1);
%The CDF is F(x)=1-((beta*a_bar)/x)^theta so we want the inverse of that to
%map a point on the CDF to a productivity draw from the distribution
productivitydraws = (beta*qual_vec).*((1-CDFvals).^(-1/theta));
